library(bugcount)
data_dir <- file.path(system.file(package = "bugcount"),"extdata")
Whitefly analysis.
Read WF count data.
WF counts per picture. Each picture must have a unique combination of
experiment + replicate + plant + pic
experiment description consists of:
exp_year + exp_cross + exp_propagation + exp_substrate
wf_file <- file.path(data_dir, "wf_consolidated.tab")
wf <- read_wf(file = wf_file, sep = "\t",header=TRUE)
Counts per Plant
wf_count <- plant_count(wf)
summary(wf_count$nymphs)
Min. 1st Qu. Median Mean 3rd Qu. Max.
0 555 1859 2506 3460 25317
What indicator should we use as resistance phenotype? Nymph counts are not normal Error is not distributed normally. Nymph counts do not meet linear model assumptions (ANOVA, MapQTL). Best fit is a negative Binomial Distribution
Probability density model fit to actual distribution
wf_count_fit <- count_fit(wf_count)
plot_count_fit(wf_count_fit, main = "Nymphs")

Reproducibility
Assume that nymphs per plant are distributed as negative binomial
Whole population.
nymphs ~ experiment GLM, posthoc and common letter difference
by_experiment <- fit_nb_glm(nymphs ~ experiment, wf_count)
plot_fit_nb_glm(wf_count,by_experiment, type = "violin", cld = FALSE)

add infestation regime as factor high,low: inferred from nyphs ~ experiment GLM
wf_count$infestation <- "low"
wf_count[wf_count$exp_year > 2016,]$infestation <- "high"
wf_count$infestation <- factor(wf_count$infestation, levels =c("low","high"))
Correlation between experiments
cor_title <- paste("Geometric Means Correlation","Log Scale", sep = "\n")
for (cross in levels(wf_count$exp_cross)) {
wf_by_exp <- wf_wide(clone ~ experiment,
wf_count[wf_count$exp_cross == cross,])
plot_wide_cor(wf_by_exp, main = paste(cross,cor_title))
}


Checks.
wf_check <- wf_count[grepl("internal_check", wf_count$group) &
wf_count$clone != "Secundina",]
wf_check <- remove_levels(wf_check)
checks_count_fit <- count_fit(wf_check)
plot_count_fit(checks_count_fit, main = "Nymphs")

Correlation between experiments
check_by_exp <- wf_wide(clone ~ experiment, wf_check,
fun = function(x) log10(geometric_mean(x)) )
plot_check_cor(check_by_exp, main = paste("Checks", cor_title))

# nymphs ~ clone GLM posthoc and common letter difference
fit <- fit_nb_glm(nymphs ~ clone, wf_check)
plot_fit_nb_glm(wf_check,fit,type = "violin", xmax = 10000)


Nymph counts by Cross

LS0tCnRpdGxlOiAiQnVnIENvdW50IEFuYWx5c2lzIgpvdXRwdXQ6CiAgaHRtbF9ub3RlYm9vazoKICAgIGhpZ2hsaWdodDogdGFuZ28KICAgIG51bWJlcl9zZWN0aW9uczogeWVzCiAgICB0aGVtZTogc3BhY2VsYWIKICAgIHRvYzogeWVzCiAgICB0b2NfZmxvYXQ6IHllcwotLS0KCmBgYHtyfQpsaWJyYXJ5KGJ1Z2NvdW50KQpkYXRhX2RpciA8LSBmaWxlLnBhdGgoc3lzdGVtLmZpbGUocGFja2FnZSA9ICJidWdjb3VudCIpLCJleHRkYXRhIikKCmBgYAoKIyBXaGl0ZWZseSBhbmFseXNpcy4KIyMgUmVhZCBXRiBjb3VudCBkYXRhLiAKV0YgY291bnRzIHBlciBwaWN0dXJlLiBFYWNoIHBpY3R1cmUgbXVzdCBoYXZlIGEgdW5pcXVlIGNvbWJpbmF0aW9uIG9mICAKZXhwZXJpbWVudCArIHJlcGxpY2F0ZSArIHBsYW50ICsgcGljICAKZXhwZXJpbWVudCBkZXNjcmlwdGlvbiBjb25zaXN0cyBvZjogIApleHBfeWVhciArIGV4cF9jcm9zcyArIGV4cF9wcm9wYWdhdGlvbiArIGV4cF9zdWJzdHJhdGUgIAoKCmBgYHtyfQp3Zl9maWxlIDwtIGZpbGUucGF0aChkYXRhX2RpciwgIndmX2NvbnNvbGlkYXRlZC50YWIiKQp3ZiA8LSByZWFkX3dmKGZpbGUgPSB3Zl9maWxlLCBzZXAgPSAiXHQiLGhlYWRlcj1UUlVFKQpgYGAKIyMgQ291bnRzIHBlciBQbGFudCAgCgpgYGB7cn0KCndmX2NvdW50IDwtIHBsYW50X2NvdW50KHdmKQpzdW1tYXJ5KHdmX2NvdW50JG55bXBocykKCmBgYAoKV2hhdCBpbmRpY2F0b3Igc2hvdWxkIHdlIHVzZSBhcyByZXNpc3RhbmNlIHBoZW5vdHlwZT8KTnltcGggY291bnRzIGFyZSBub3Qgbm9ybWFsIApFcnJvciBpcyBub3QgZGlzdHJpYnV0ZWQgbm9ybWFsbHkuCk55bXBoIGNvdW50cyBkbyBub3QgbWVldCBsaW5lYXIgbW9kZWwgYXNzdW1wdGlvbnMgKEFOT1ZBLCBNYXBRVEwpLgpCZXN0IGZpdCBpcyBhIG5lZ2F0aXZlIEJpbm9taWFsIERpc3RyaWJ1dGlvbgoKIyMgUHJvYmFiaWxpdHkgZGVuc2l0eSBtb2RlbCBmaXQgdG8gYWN0dWFsIGRpc3RyaWJ1dGlvbgpgYGB7ciwgd2FybmluZz1GQUxTRX0Kd2ZfY291bnRfZml0IDwtICBjb3VudF9maXQod2ZfY291bnQpCgpwbG90X2NvdW50X2ZpdCh3Zl9jb3VudF9maXQsIG1haW4gPSAiTnltcGhzIikKYGBgCgojIyBSZXByb2R1Y2liaWxpdHkgCgpBc3N1bWUgdGhhdCBueW1waHMgcGVyIHBsYW50IGFyZSBkaXN0cmlidXRlZCBhcyBuZWdhdGl2ZSBiaW5vbWlhbAoKIyMjIFdob2xlIHBvcHVsYXRpb24uCm55bXBocyB+IGV4cGVyaW1lbnQgR0xNLCBwb3N0aG9jIGFuZCBjb21tb24gbGV0dGVyIGRpZmZlcmVuY2UgCgpgYGB7cn0KYnlfZXhwZXJpbWVudCA8LSBmaXRfbmJfZ2xtKG55bXBocyB+IGV4cGVyaW1lbnQsIHdmX2NvdW50KQoKYGBgCgpgYGB7cn0KcGxvdF9maXRfbmJfZ2xtKHdmX2NvdW50LGJ5X2V4cGVyaW1lbnQsIHR5cGUgPSAidmlvbGluIiwgY2xkID0gRkFMU0UpCmBgYAphZGQgaW5mZXN0YXRpb24gcmVnaW1lIGFzIGZhY3RvciAKaGlnaCxsb3c6IGluZmVycmVkIGZyb20gbnlwaHMgfiBleHBlcmltZW50IEdMTQoKYGBge3J9CndmX2NvdW50JGluZmVzdGF0aW9uIDwtICJsb3ciCndmX2NvdW50W3dmX2NvdW50JGV4cF95ZWFyID4gMjAxNixdJGluZmVzdGF0aW9uIDwtICJoaWdoIgp3Zl9jb3VudCRpbmZlc3RhdGlvbiA8LSBmYWN0b3Iod2ZfY291bnQkaW5mZXN0YXRpb24sIGxldmVscyA9YygibG93IiwiaGlnaCIpKQpgYGAKCiMjIENvcnJlbGF0aW9uICBiZXR3ZWVuIGV4cGVyaW1lbnRzIAoKYGBge3IsIGZpZy5hc3A9MX0KCmNvcl90aXRsZSA8LSBwYXN0ZSgiR2VvbWV0cmljIE1lYW5zIENvcnJlbGF0aW9uIiwiTG9nIFNjYWxlIiwgc2VwID0gIlxuIikKCmZvciAoY3Jvc3MgaW4gbGV2ZWxzKHdmX2NvdW50JGV4cF9jcm9zcykpIHsKICAKICB3Zl9ieV9leHAgPC0gd2Zfd2lkZShjbG9uZSB+IGV4cGVyaW1lbnQsCiAgICAgICAgICAgICAgICAgICAgICAgd2ZfY291bnRbd2ZfY291bnQkZXhwX2Nyb3NzID09IGNyb3NzLF0pCiAgCiAgcGxvdF93aWRlX2Nvcih3Zl9ieV9leHAsIG1haW4gPSBwYXN0ZShjcm9zcyxjb3JfdGl0bGUpKQogIAp9CmBgYAoKIyMgQ2hlY2tzLgpgYGB7ciwgd2FybmluZz1GQUxTRX0Kd2ZfY2hlY2sgPC0gd2ZfY291bnRbZ3JlcGwoImludGVybmFsX2NoZWNrIiwgd2ZfY291bnQkZ3JvdXApICYKICAgICAgICAgICAgICAgICAgICAgICB3Zl9jb3VudCRjbG9uZSAhPSAiU2VjdW5kaW5hIixdCgp3Zl9jaGVjayA8LSByZW1vdmVfbGV2ZWxzKHdmX2NoZWNrKQoKY2hlY2tzX2NvdW50X2ZpdCA8LSAgY291bnRfZml0KHdmX2NoZWNrKQoKcGxvdF9jb3VudF9maXQoY2hlY2tzX2NvdW50X2ZpdCwgbWFpbiA9ICJOeW1waHMiKQpgYGAKCgojIyMgQ29ycmVsYXRpb24gIGJldHdlZW4gZXhwZXJpbWVudHMKCmBgYHtyLCBmaWcuYXNwPTF9CmNoZWNrX2J5X2V4cCA8LSB3Zl93aWRlKGNsb25lIH4gZXhwZXJpbWVudCwgd2ZfY2hlY2ssIAogICAgICAgICAgICAgICAgICAgICAgICBmdW4gPSBmdW5jdGlvbih4KSBsb2cxMChnZW9tZXRyaWNfbWVhbih4KSkgKQpwbG90X2NoZWNrX2NvcihjaGVja19ieV9leHAsIG1haW4gPSBwYXN0ZSgiQ2hlY2tzIiwgY29yX3RpdGxlKSkKYGBgCgoKYGBge3J9CiMgbnltcGhzIH4gY2xvbmUgR0xNIHBvc3Rob2MgYW5kIGNvbW1vbiBsZXR0ZXIgZGlmZmVyZW5jZQpmaXQgPC0gZml0X25iX2dsbShueW1waHMgfiBjbG9uZSwgd2ZfY2hlY2spCnBsb3RfZml0X25iX2dsbSh3Zl9jaGVjayxmaXQsdHlwZSA9ICJ2aW9saW4iLCB4bWF4ID0gMTAwMDApCmBgYAoKYGBge3IsIGZpZy5hc3A9MC41fQpmaXQgPC0gTUFTUzo6Z2xtLm5iKG55bXBocyB+ICBpbmZlc3RhdGlvbiAqIGNsb25lLCBkYXRhID0gd2ZfY2hlY2spCnBvc3Rob2MgPC0gbHNtZWFuczo6Y2xkKGxzbWVhbnM6OmxzbWVhbnMoZml0LCB+IGluZmVzdGF0aW9uICogY2xvbmUsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgYWRqdXN0ID0gInR1Y2tleSIpLCB0eXBlID0gInJlc3BvbnNlIikKIyBQbG90IHBvc3QgaG9jCgojIyMgT3JkZXIgdGhlIGxldmVscyBmb3IgcHJpbnRpbmcKY2xvbmVfb3JkZXIgPC0gd2ZfYWdncmVnYXRlKHdmX2NoZWNrLCBieV9zdGF0ID0gIm1lYW4iKVsiY2xvbmUiXQoKcG9zdGhvYyRjbG9uZSA8LSBmYWN0b3IocG9zdGhvYyRjbG9uZSwgbGV2ZWxzID0gY2xvbmVfb3JkZXIkY2xvbmUpCgpwb3N0aG9jJGluZmVzdGF0aW9uID0gZmFjdG9yKHBvc3Rob2MkaW5mZXN0YXRpb24sCiAgICAgICAgICAgICAgICAgICBsZXZlbHM9YygibG93IiwgImhpZ2giKSkKCiMjIyAgUmVtb3ZlIHNwYWNlcyBpbiAuZ3JvdXAgIAoKcG9zdGhvYyQuZ3JvdXA9Z3N1YigiICIsICIiLCBwb3N0aG9jJC5ncm91cCkKCiMgcXVhcnR6KCkKcGxvdF9nZ19pbmZlc3RhdGlvblhjbG9uZShwb3N0aG9jKQpgYGAKCgoKIyBOeW1waCBjb3VudHMgYnkgQ3Jvc3MgCmBgYHtyLCBmaWcuYXNwPTF9CmV4cF9sZXZlbHMgPC0gbGV2ZWxzKGFzLmZhY3Rvcih3Zl9jb3VudCRleHBlcmltZW50KSkKZXhwX2dycCA8LSBsaXN0KENNODk5NiA9IGxpc3QoMSwyLDQsNSksCiAgICAgICAgICAgICAgICBHTTg1ODYgPSBsaXN0KDMsNiw3KSkKcGxvdF9saXN0IDwtIGxpc3QoKQpmb3IgKGNyb3NzIGluIGxldmVscyh3Zl9jb3VudCRleHBfY3Jvc3MpKSB7CiAgZm9yIChpZHggaW4gMTpsZW5ndGgoZXhwX2dycFtbY3Jvc3NdXSkpIHsKICBleHBfYWxsb3dlZCA8LSBleHBfbGV2ZWxzW2V4cF9ncnBbW2Nyb3NzXV1bW2lkeF1dXQogIAogIHdmX2J5X2Nyb3NzIDwtIHdmX2NvdW50W3dmX2NvdW50JGV4cF9jcm9zcyA9PSBjcm9zcyAmCiAgICAgICAgICAgICAgICAgICAgICAgICAgICB3Zl9jb3VudCRleHBlcmltZW50ICVpbiUgZXhwX2FsbG93ZWQsXQogIGV4cF9jb3VudCA8LSBhZ2dyZWdhdGUoZXhwZXJpbWVudCB+IGNsb25lICx3Zl9ieV9jcm9zcywKICAgICAgICAgICAgICAgICAgICAgICAgIEZVTiA9IGZ1bmN0aW9uKHgpIGxlbmd0aCh1bmlxdWUoeCkpKQogIHNlbGVjdGVkX2Nsb25lcyA8LSB3aXRoKGV4cF9jb3VudCwgewogICAgY2xvbmVbZXhwZXJpbWVudCA9PSBsZW5ndGgoZXhwX2FsbG93ZWQpXQogIH0pCiAgCiAgd2ZfYnlfY3Jvc3MgPC0gd2ZfYnlfY3Jvc3NbIHdmX2J5X2Nyb3NzJGNsb25lICVpbiUgc2VsZWN0ZWRfY2xvbmVzLF0KICB3Zl9ieV9jcm9zcyA8LSAod2ZfYnlfY3Jvc3MpCiAgCiAgd2ZfYnlfY3Jvc3MkbnltcGhzIDwtIHdmX2J5X2Nyb3NzJG55bXBocysxCiAgIyBHTE0gKioqKioqKioqKioqKioqKioqKioqKioqKioqKioqKioqKioqKioqKioqKioqKioqKioqKioqKioqKioqICMKICAjIEFzc3VtaW5nIE5lZ2F0aXZlIEJpbm9taWFsIGRpc3RyaWJ1dGlvbgogICMKICAKICBmaXRfbmIgPC0gTUFTUzo6Z2xtLm5iKG55bXBocyB+ICBjbG9uZSwgZGF0YSA9IHdmX2J5X2Nyb3NzKQogIEFJQyhmaXRfbmIpCgogIHBvc3Rob2MgPC0gbHNtZWFuczo6Y2xkKAogIGxzbWVhbnM6OmxzbWVhbnMoCiAgICAgZml0X25iLCB+IGNsb25lLCBhZGp1c3QgPSAidHVja2V5IiksdHlwZSA9ICJyZXNwb25zZSIKICApCiAgCiAgIyBQbG90IHBvc3QgaG9jCiAgYnJlYWtzIDwtIHBvc3Rob2MkY2xvbmVbIWdyZXBsKCJHTXxDTSIsIHBvc3Rob2MkY2xvbmUsIHBlcmwgPSBUUlVFKV0KICAKICBwbG90X2xpc3RbW3Bhc3RlKGV4cF9hbGxvd2VkKV1dIDwtIHBsb3RfY3Jvc3NfbmJfZml0KHBvc3Rob2MsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBleHBfYWxsb3dlZCwKICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGJyZWFrcyA9IGJyZWFrcykKICAKICB9Cn0KZ2dwbG90Mjo6dGhlbWVfc2V0KGdncGxvdDI6OnRoZW1lX2dyYXkoYmFzZV9zaXplID0gMTApKQpvcmRlcmVkX2xpc3QgPC0gYyhwbG90X2xpc3RbMV0sbGlzdChlbXB0eSA9IE5VTEwpLHBsb3RfbGlzdFtjKDMsNywyLDUsNCw2KV0pCm5hbWVzKHBsb3RfbGlzdCkKbmFtZXMob3JkZXJlZF9saXN0KQpjb3dwbG90OjpwbG90X2dyaWQocGxvdGxpc3QgPSBvcmRlcmVkX2xpc3QsCiAgICAgICAgICBuY29sID0gMiwgbnJvdyA9IDQsIGFsaWduID0gInYiKQoKYGBgCg==